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Abstract: Soil erosion in the Three-River Headwaters Region (TRHR) of the Qinghai-Tibet Plateau in 
China has a significant impact on local economic development and ecological environment. Vegetation 
and precipitation are considered to be the main factors for the variation in soil erosion. However, it is a 
big challenge to analyze the impacts of precipitation and vegetation respectively as well as their combined 
effects on soil erosion from the pixel scale. To assess the influences of vegetation and precipitation on the 
vatiation of soil erosion from 2005 to 2015, we employed the Revised Universal Soil Loss Equation 
(RUSLE) model to evaluate soil erosion in the TRHR, and then developed a method using the 
Logarithmic Mean Divisia Index model (LMDI) which can exponentially decompose the influencing 
factors, to calculate the contribution values of the vegetation cover factor (C factor) and the rainfall 
erosivity factor (R factor) to the variation of soil erosion from the pixel scale. In general, soil erosion in 
the TRHR was alleviated from 2005 to 2015, of which about 54.95% of the area where soil erosion 
decreased was caused by the combined effects of the C factor and the R factor, and 41.31% was caused by 
the change in the R factor. There were relatively few areas with increased soil erosion modulus, of which 
64.10% of the area where soil erosion increased was caused by the change in the C factor, and 23.88% was 
caused by the combined effects of the C factor and the R factor. Therefore, the combined effects of the 
C factor and the R factor were regarded as the main driving force for the decrease of soil erosion, while 
the C factor was the dominant factor for the increase of soil erosion. The area with decreased soil erosion 
caused by the C factor (12.1010? km?) was larger than the area with increased soil erosion caused by the 
C factor (8.30X10° km), which indicated that vegetation had a positive effect on soil erosion. This study 
generally put forward a new method for quantitative assessment of the impacts of the influencing factors 
on soil erosion, and also provided a scientific basis for the regional control of soil erosion. 
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1 Introduction 


Soil erosion is a worldwide problem with varying degrees of impacts on the life of human beings 
and has aroused widespread attention (Thomas et al., 2018). The Three-River Headwaters Region 
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(TRHR) of the Qinghai-Tibet Plateau in China is a crucial barrier for ecological protection with 
some important ecosystem services, such as water and soil conservation and biodiversity 
protection (Jiang et al., 2016; Shao et al., 2017). In recent years, global warming and human 
activities have affected the region's ecological environment to some extent, which subsequently 
led to the rise of snow lines, the reduction of runoff, the increasing frequency of river closure, and 
the changes in vegetation coverage and biomass (Cook et al., 2016). Furthermore, located in the 
hinterland of the Qinghai-Tibet Plateau, the TRHR can hardly be recovered once it was degraded, 
as this region is a typical area with a fragile ecological environment characterized by poor 
climatic conditions, severe soil erosion, harsh natural conditions, and frequent natural disasters. 
Before the implementation of the ecological conservation and restoration project, the alpine 
meadows in the TRHR had shown a complete degradation trend (Zheng et al., 2019). 
Consequently, soil erosion in the TRHR is severe, with the total potential amount of 1.12x10° t/a 
(Liu et al., 2005). The data illustrated that the area of soil erosion in this region in 2010 was as 
high as 114.80x10? km’, accounting for more than 30% of the total area, which seriously affected 
the ecological security of the region. Soil erosion is affected by vegetation, precipitation, soil 
property, topography, and human activity (Thomas et al., 2018). The change in each factor can 
cause varying degrees of impacts on soil erosion. Generally speaking, soil property and 
topography are relatively stable with little changes. Therefore, the variation of soil erosion is 
mainly driven by precipitation and vegetation. Due to their interactions, the relationships among 
precipitation, vegetation, and soil erosion are uncertain and complex. Hence, it is significant to 
understand the temporal and spatial changes of soil erosion and the impacts of the influencing 
factors on soil erosion in the TRHR. 

The formation mechanism and the identification of the influencing factors of soil erosion are 
the cores and frontier issues of current researches. However, much more efforts should be put into 
the quantitative attribution study of the combined effects of multiple factors (Nearing et al., 2004; 
Zuazo and Pleguezuelo, 2009; Panagos et al., 2015a; Ochoa et al., 2016; Guerra et al., 2017). So 
far, extensive researches on soil erosion and its influencing factors have been carried out by both 
domestic and foreign scholars. For examples, Ganasri et al. (2016) studied soil erosion in 
Nethravathi Basin using remote sensing data based on the Revised Universal Soil Loss Equation 
(RUSLE) model; Guerra et al. (2016) gave a comprehensive introduction and a review of the 
current trends of ecosystem services provision by assessing soil erosion prevention in 
Mediterranean areas. Studies on the effects of rainfall patterns on runoff and soil erosion in field 
plots were conducted by Mohamadi et al. (2015). Garcia-Ruiz et al. (2010) assessed the effects of 
land use on soil erosion in Spain. Panagos et al. (2015a) proposed a new European slope length 
and steepness factor to model soil erosion caused by water. According to these researches, the 
main methods used in the study of influencing factors of soil erosion are as follows: principal 
component analysis, regression analysis, and correlation coefficient method. These methods can 
determine the correlation between these factors and soil erosion from various aspects in different 
regions, and estimate the impacts of different factors on soil erosion. However, they cannot 
quantitatively judge the specific contribution value of each factor to the variation in soil erosion 
from the pixel scale. At present, there are few pieces of research to discuss the contribution values 
of various influencing factors to the variation in soil erosion from the pixel scale. It is meaningful 
and essential to scientifically explore the variation in soil erosion affected by various factors from 
the pixel scale, which is conducive to the ecological protection and the control of soil erosion in 
local areas. 

Logarithmic Mean Divisia Index (LMDI) method is a preferred approach for quantifying the 
impact of different factors on the change of aggregates because of its profound theoretical 
foundation, adaptability, usability, and some other desirable properties in the context of 
decomposition analysis (Ang, 2004). It can be quite useful to find out the influencing factors of 
aggregate indicator, analyze the influencing degree of each factor by decomposing the aggregate 
indicator, and explain the reasons for the change in aggregate indicator (Ang, 2005, 2015). 
Compared with other decomposition models, this method provides a logarithmic mean weight 
equation without residuals, which can give a perfect decomposition without any unexplained 
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residual term appearing in the results. Also, it has the advantage of decomposing quantity 
indicator as well as intensity indicator (Ang, 2004; Goh and Ang, 2018). As an analytical 
framework for studying the characteristics and the acting mechanism of the variation, the index 
decomposition analysis method has been frequently utilized in social and economic researches in 
recent years (Wang and Ang, 2018). Because of its unique advancement and flexibility, the 
method has been successfully applied in the field of analyzing the driving factors for energy 
consumption and CO; emissions (Mousavi et al., 2017; Moutinho et al., 2018; Pourebadollahan 
et al., 2018). It mainly decomposes the total increment from different layers and quantifies the 
contribution of structural changes from different levels to the total increment, in which it has 
achieved remarkable results. In addition, the method has been successfully applied in the studies 
of water resources (Cazcarro et al., 2019; Llop, 2019). LMDI can completely decompose the 
multiple related factors at the same time. It has passed the factor variance test and the time 
variance test as well as solved the problem of zero value in data. LMDI has many advantages, 
for instance, it has no residual error and can also be very convenient to use. The empirical 
RUSLE, which can predict the average annual soil loss resulting from raindrop and runoff from 
field slopes (Renard, 1997; Panagos et al., 2015b), is currently widely accepted as a method for 
evaluating soil and water conservation functions. The equation, which is in line with the LMDI 
model, is composed of vegetation cover factor, soil erodibility factor, rainfall erosivity factor, 
support practice factor, and topographic factor. In this study, the LMDI model is first introduced 
and applied to estimate the contribution values of different factors to the variation of soil 
erosion. 

In the TRHR, the first-stage of ecological conservation and restoration project was launched in 
2005 and completed in 2013. During the period, large-scale vegetation construction and 
restoration efforts were carried out, resulting in significant changes in the overall ecological 
environment in the region. How did the spatial differentiation of soil erosion variation as a result 
of climate changes and increased human activities? How did soil erosion respond to the 
ecological project and climate changes? It is of great research value to clarify these scientific 
issues for policymaking and adjustments relating to ecological environment governance in the 
TRHR. Therefore, on the basis of the RUSLE and LMDI models, this study will quantitatively 
evaluate the dynamic changes in soil erosion and the effects of the influencing factors on soil 
erosion in the TRHR. Our purpose is to provide a new method to analyze the influencing factors 
of soil erosion using the LMDI model, as well as to put forward a reference for policy 
formulation to control soil erosion in the TRHR. 


2 Materials and methods 


2.1 Study area 


The study area (TRHR) is located in the hinterland of the Qinghai-Tibet Plateau (31°32'—36°17'N, 
89°24’-102°15’'E; Fig. 1), China. The TRHR is the birthplace of the Yangtze River, the Yellow 
River, and the Lancang River. It is not only an important ecological functional zone for water 
conservation in China but also one of the most sensitive and vulnerable regions of the ecosystem 
(Nie et al., 2018). The study area is vast, covering approximately 350.6010? km? with complex 
terrains. It is dominated by mountainous landforms with an average elevation of about 4592.87 m. 
The region is characterized by a typical plateau continental climate, with alternating hot and cold 
seasons, distinct wet and dry seasons, long sunshine hours, and strong radiation. The annual mean 
temperature is between —5.6°C and 7.8°C, and the mean annual precipitation is between 262.2 
and 772.8 mm. The spatial distribution of precipitation is uneven, which gradually decreases from 
the southeast to the northwest, with significant regional differentiation. The soil in the study area 
is barren, and the main soil type is alpine meadow soil. The permafrost is extremely developed in 
this region (Guo et al., 2017). In the study area, there are various types of natural vegetation 
types, such as coniferous forests, shrubs, alpine meadows, alpine grasslands, and alpine sparse 
vegetation, which are distributed successively from the southeast to the northwest. 
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Fig. 1 Location of the Three-River Headwaters Region (TRHR) in the Qinghai-Tibet Plateau, China (a) and 
overview of the TRHR (b). DEM, digital elevation model. 


2.2 Data and processing 
2.2.1 Normalized Difference Vegetation Index (NDVI) data 


The NDVI data used in the study were collected from the Resource and Environment Data Cloud 
Platform (http://www.resdc.cn/), including SPOT/VEGETATION NDVI datasets in 2005 and 
2015. The spatial resolution is 1 km and the time resolution is 10 d. The digital images were 
processed by the Maximum Value Composite (MVC) method to extract the maximum value of 
NDVI in each year. The MVC algorithm is the most common compositing criterion used to 
generate NDVI composite data (Gu et al., 2009), which can eliminate noise from NDVI series as 
well as errors due to sensor-related artefacts, resulting in more reliable NDVI values (Beck et al., 
2007). After image mosaic and clip, we obtained the NDVI composite data of the study area in 
2005 and 2015. 

2.2.2 Meteorological observation data 


Meteorological data in 2005 and 2015 were provided by the National Meteorological Center of 
China (http://data.cma.cn/) from 174 meteorological stations, some of which were in the study 
area and others were located in the vicinity of the study area (Fig. 1a). The data outliers were 
removed before calculating the rainfall erosivity. Studies have shown that traditional interpolation 
approaches such as Inverse Distance Weighting (IDW) and Kriging methods performed not very 
well on spatial interpolation of meteorological data, especially in some highly heterogeneous 
areas such as the Qinghai-Tibet Plateau (Tan et al., 2016; Sha et al., 2017). Thus, we used the 
Australian Smooth Spline Function Interpolation Tool (ANUSPLIN), a recognized professional 
interpolation method (Liu et al., 2008) to interpolate meteorological data in this study. 
ANUSPLIN is a multi-variable climate interpolation tool for smoothing spline function and it 
takes elevation as the covariate in interpolation (Hutchinson and Xu, 2004; Plouffe et al., 2015), 
which has been demonstrated to have a higher accuracy in the interpolation of meteorological 
data. In this study, the interpolation of rainfall was carried out by using ANUSPLIN, and the 
interpolated results of 18 meteorological stations which are completely distributed within the 
TRHR were compared with the observed values (Fig. 2). The correlation coefficient (R*) of both 
years (2005 and 2015) was greater than or equal to 0.90, indicating that the interpolation results 
were good. The spatial distribution of total annual rainfall in 2005 and 2015 in the TRHR are 
shown in Figure S1 in Supplementary material. 

2.2.3 Soil and digital elevation model (DEM) data 


Soil attribute data were obtained from the soil database of China (Shi et al., 2004) with a scale of 
1:1,000,000. The database includes the spatial distribution of soil types, the physical properties of 
various soil types, the structure of soil, and the percentage contents of sand, silt, clay and organic 
matter. Soil types in the study area were screened and matched with the corresponding soil 
properties in the database to calculate the erodibility factor. The DEM data were derived from the 
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Fig. 2 Verification of the annual rainfall interpolation results in 2005 (a) and 2015 (b) 


Shuttle Radar Topography Mission images in Geo-TIFF format provided by the Geospatial Data 
Cloud (http://www.gscloud.cn/) with a resolution of 90 m and then were resampled to a resolution 
of 1 km for this study. The DEM outliers were removed using the ArcGIS platform, and the 
spatial distribution data of the topographic factor for the study area were extracted from the 
processed DEM. 

2.2.4 Land use data 

The land use data of 2005 and 2015 were collected from the Resource and Environment Data 
Cloud Platform (http://www.resdc.cn/). The classification of land use types is based on the 
second-level classification standard of the Chinese Academy of Sciences. According to the natural 
attributes of land resources, there are 25 types of land use, such as paddy fields, drylands and 
lakes. The land use map of the study area was used to assign support practice (P) factor values. 
2.2.5 Statistic data 

The soil and water conservation statistics in 2015 were obtained from Qinghai Provincial Water 
Resources Department (http://slt.qinghai.gov.cn/), which recorded the amount of soil erosion in 
the TRHR in 2015. The data were used to verify the accuracy of simulated soil erosion calculated 
by the RUSLE model in this study. 

It should be noted that we uniformly projected the data using the Krasovsky_1940_Albers 
coordinate system in order to ensure that data used for this study had a good spatial consistency. 
As for the resolution, in consideration of the 1 km spatial resolution of the remote sensing data 
and the land use data, we resampled all data to a common spatial resolution of 1 kmx1 km. 


2.3 Methods 


The database for the assessment of soil erosion and the analysis of influencing factors was 
constructed by collecting data such as vegetation, meteorology, and topography. The workflow and 
the methods for the study are presented in Figure 3. First, we evaluated soil erosion in 2005 and 
2015 based on the RUSLE model in the study area. Specifically, we used the statistics from 
Qinghai Provincial Water Resources Department (http://slt.qinghai.gov.cn/) and the results from 
previous researches to verify the simulated results. Then, we analyzed the spatial distribution and 
variation characteristics of soil erosion during the study period. Next, on the basis of the 
assessment results of soil erosion, we designed a method based on the LMDI model to 
quantitatively analyze the contribution values of the vegetation cover (C) factor and rainfall 
erosivity (R) factor to the variation in soil erosion from the pixel scale. Finally, we analyzed the 
spatial distribution characteristics of the relative impacts of each influencing factor on soil erosion. 
2.3.1 Revised universal soil loss equation (RUSLE) model 

Wischmeier and Smith (1978) established the universal soil loss equation (USLE) based on a 
large number of regional observations and artificially simulated rainfall experiments; the US 
Department of Agriculture's Agricultural Research Service (USDA-ARS) put forward an 
advanced assessment model for soil erosion, namely the Revised Universal Soil Loss Equation 
(RUSLE) model, based on the original USLE model. RUSLE is one of the most commonly used 
soil erosion equations. The model estimates soil erosion caused by five factors that represent the 
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Fig. 3 Workflow and methods of the study. DEM, digital elevation model; RUSLE, revised universal soil loss 
equation; LMDI, logarithmic mean Divisia index; A, the amount of soil erosion; C factor, vegetation cover factor; 
R factor, rainfall erosivity factor; LS factor, topographic factor; K factor, soil erodibility factor; P factor, support 
practice factor. 


effects of vegetation, precipitation, soil properties, topography, and land use on soil erosion. The 
average annual soil loss was obtained in RUSLE by multiplying the factors manifested in the 
following equation: 

A=CxRXLSXxK*P, (1) 
where A is the amount of soil erosion (t/(hm?.a)); C stands for the vegetation cover factor 
(dimensionless), indicating the influence of vegetation coverage on soil erosion, which can be 
calculated by the method of Cai et al. (2000); and R is the rainfall erosivity factor 
(MJ-mm/(hm?-h-a)), indicating the effects of precipitation on soil erosion. The R factor was 
computed using the equation developed by Wischmeier and Smith (1978). LS is the topographic 
factor (dimensionless), indicating the effects of steepness and slope length on soil erosion. The LS 
factor (dimensionless) was obtained by the method proposed by Fu et al. (2015), which is suitable 
for China. K refers to the soil erodibility factor (t-hm?-h/(hm?-MJ-mm)), which represents the 
effects of soil types on soil erosion and was calculated by the EPIC (erosion-productivity impact 
calculator) model (Sharpley and Williams, 1990). P stands for the support practice factor 
(dimensionless), which is defined as the ratio of soil loss with a specific support practice to the 
corresponding loss with up slope and down slope cultivation (Dissanayake et al., 2019). The value 
of P factor ranges from 0 to 1, in which O represents high-quality conservation practices and 1 
indicates poor preservation conservation practices. According to the previous study of soil erosion 
in the Qinghai-Tibet Plateau (Kang et al., 2018), the P factor value of agricultural lands in the 
study area is 0.15. There is no erosion in the water, wetlands, bare rocks, snow, and ice, thus the P 
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factor value equals to 0; and, the P value equals to 1 in the remaining land use types since there is 
no any water conservation measure. The spatial distribution maps of each factor calculated by the 
methods above are shown in Figure 4. The specific calculation methods of the RUSLE factors are 
presented in Section 2 in Supplementary material. 
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Fig.4 Layers of the RUSEL factors in the TRHR. (a), distribution of C factor values in 2005; (b), distribution of 
C factor values in 2015; (c), distribution of K factor values; (d), distribution of LS factor values; (e), distribution 
of P factor values in 2005; (f), distribution of P factor values in 2015; (g), distribution of R factor values in 2005; 
(hb), distribution of R factor values in 2015. 


2.3.2 LMDI model 


As described in the foregoing, the LMDI is the preferred method in decomposition analysis. The 
general structure of the decomposition of the LMDI model is as follows (Ang, 2005). 
The target variable is assumed as V. The changes in V are related to n index factors over time. 
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Each factor is linked with n index variables (x1, X2, X3, ..., Xn). The general index decomposition 
analysis is shown in Equation 2. 


VOS (2) 


where the subscript i represents a sub-category level of the target variable. 
In addictive decomposition, the difference of the target variable is caused by respective factors 
(AVx, AVx2, ..., AVxn). 


AV, 


ot =V" -V° = AVa + AV,, ++ AV, (3) 
where AVio represents the total changes in the target variable; V? stands for the target variable in 
the base period 0; V” stands for the target variable in the planned period T; and AV x (k=1, 2, 3, ..., 
n) denotes the difference of the target variable associated with the respective factors. 

In the multiplicative decomposition, the ratio caused by respective factors (Du, Dx2, ..., Dxn) 


was decomposed as follows: 

Dep =V LV? = Dy X Dy XX Digs (4) 
where D,.: represents the ratio of the target variable between the planned period T and the base 
period 0. Dx (k=1, 2, 3, ..., n) stands for the ratio of the change in the target variable associated 
with the respective factors (X1, X2, X3, ..., Xn). 

Using the LMDI approach proposed by Ang et al. (1998), the formula used to calculate the 
contribution value of each factor to the change of the target variable was derived, as shown in 
Equation 5. In the multiplicative decomposition, the contribution values of each factor to the 
variation of the target variable are shown in Equations 6 and 7. 
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where AV,x and D,x are the effects of the k factor to the target variable V in the addictive and 
multiplicative decomposition, respectively; xi; is the value of the k® factor in the base period 0; 
Xi: is the value of the k" factor in the planned period T; a and b are two positive numbers; and 
L(a,b) is the function representing the logarithmic average of a and b. 

Comparing the LMDI approach and RUSLE model, we can find that the structure of the 
RUSLE conforms to that of the LMDI. The amount of soil erosion is the target variable and the 
index factors are C, R, LS, K, and P, which means that n is equal to 5 and the variables are xı, x2, 
X3, X4, and Xs, respectively. 

Since each factor in the RUSLE model is a spatially distributed raster image dataset, the 
contribution values of each factor to soil erosion modulus for each pixel can be calculated using 
the LMDI approach. Additive decomposition and multiplicative decomposition were equally valid 
and had the equivalent interpretive power (Ang, 2005, 2015). Compared with the multiplicative 
case, the decomposition results of the additive case are given as physical units instead of indices 
so that the addictive decomposition is easier to be interpreted and utilized (Zhang et al., 2018). 
2.3.3 Contributions of the influencing factors to the variation of soil erosion 
It can be seen from the RUSLE equations that the main factors influencing soil erosion are 
topography, soil, vegetation, precipitation, and support practices. It should be noted that we did not 
evaluate the contribution values of the LS factor, the K factor, and the P factor to the variation of 
soil erosion due to the following reasons. First, the LS factor barely changed during the study 
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period. Second, on top of the little changes in the K factor mentioned above, the soil dataset used 
to calculate the K factor value was from the Second National Soil Survey, which is the only one 
dataset. Thus, we did not consider the change of the K factor due to data limitation. Third, for the P 
factor, it was mainly determined by the land use types. Because there was little variation in land 
use types from 2005 to 2015 (area percentage of cropland increased by 0.0009%; area percentage 
of water, wetland, bare land, snow, and ice increased by 0.0083%; and area percentage of the other 
land use types decreased by 0.0092%), we also did not explore the impact of the P factor. 

We introduced the LMDI model to devise a method to quantitatively analyze the contribution 
values of the C factor and the R factor to the variation of soil erosion during the period when the 
ecological conservation and restoration project was implemented in the TRHR. 

According to the additive decomposition of the LMDI approach and the RUSLE model, we 
calculated the contribution values of the C factor to the variation of soil erosion from the pixel 
scale. The equation is described as follows (Eq. 8): 

(Cals 


2015 2005 In Cy; An Ae” In Gre 
AA, = SLA A )In Cae -r ApS In A cee a (8) 
where AAc is the contribution value of the C factor to the variation of soil erosion (t/hm-a)); Ai” 
is the soil erosion modulus in 2015 (t/(hm?.a)); A;™ is the soil erosion modulus in 2005 
(t/(hm?-a)); Cx; is the vegetation cover factor in 2015; and C;; is the vegetation cover factor in 
2005. The value of AAc higher than zero indicates that the C factor has aggravated soil erosion, 
while the value of AAc lower than zero indicates that the C factor has alleviated soil erosion. 
Since there are no subsets under each factor in the RUSLE model, the value of n equals to 1 in the 
equation. 
The contribution of the R factor to the variation of soil erosion can be determined by the 
following equation: 


2015 2005 ee A= LEE Re 
AA, = SLA A )in Rams -X In Aaa m In Aa In Ri. a (9) 


where AAR is the contribution value of the R factor to the variation of soil erosion (t/(hm?-a)); Ri;"” 
is the rainfall erosivity factor in 2015; and R” is the rainfall erosivity factor in 2005. The value 
of AAp higher than zero indicates that the R factor has positive effects on the increase of soil 
erosion, while the value of AAr lower than zero indicates that the R factor has negative effects on 
the increase of soil erosion. 

The spatial distribution maps of the contributions of the C factor and the R factor to the 
variation of soil erosion were obtained through quantitative calculation from the pixel scale using 
the methods mentioned above. The dominant factors that influence the increase or decrease of soil 
erosion can be identified by the spatial distribution maps. In order to further grasp the relative 
effects of the C factor and the R factor on the variation of soil erosion, we combined the values of 
AA, AAc, and AAg, resulting in six schemes in total (Table 1). 


3 Results 
3.1 Validation of the RUSLE model results 


To verify the accuracy of the simulated results of the RUSLE model, we compared the results 
with statistics from Qinghai Provincial Water Resources Department (http://slt.qinghai.gov.cn/) in 
2015. It manifested that soil erosion area estimated by the RUSLE model was 12.13x10* km?, and 
the statistic showed that the area of soil erosion was 13.07x10* km?. The ratio between the 
simulated value and the measured value was 0.93, indicating that the simulated results estimated 
by the RUSLE model was consistent with the actual statistics. Besides, we also verified the 
estimated results using the spatial distribution map of soil erosion within the study area in 2015 
obtained from the Ministry of Water Resources of the People's Republic of China 
(http://www.mwr.gov.cn/sj/tjgb/zgstbcgb/201612/t20161229_783346.html), and good consistency 
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Table 1 Criteria for assessing the relative effects of the C factor and R factor on the variation of soil erosion 


AA AAc AAR Description 

>0 >0 Increased soil erosion due to the C factor and the R factor (IECR) 
>0 >0 <0 Increased soil erosion due to the C factor (IEC) 

<0 >0 Increased soil erosion due to the R factor (IER) 

<0 <0 Decreased soil erosion due to the C factor and the R factor (DECR) 
<0 <0 >0 Decreased soil erosion due to the C factor (DEC) 

>0 <0 Decreased soil erosion due to the R factor (DER) 


Note: C factor, vegetation cover factor; R factor, rainfall erosivity factor; AA, the variation of soil erosion; AAc, contribution value of 
the C factor to the variation of soil erosion; AAp, contribution value of the R factor to the variation of soil erosion. 


was shown in the same region. The slight differences between the simulated results estimated by 
the RUSLE and the measured results or the results from previous studies in the region can be 
attributed to different data and processing methods used in the evaluation. For example, the 
methods for interpolation of meteorological data were varying. Overall, the results proved that 
soil erosion values estimated by the RUSLE model had a quite high accuracy. 


3.2 Spatial distribution and variation of soil erosion 


The spatial distribution maps of soil erosion in the study area in 2005 and 2015 were obtained 
from the RUSLE model (Fig. 5). Soil erosion modulus was classified into five classes: <100, 
100-200, 200-300, 300-400 and >400 t/(hm?.a). The results manifested that the average soil 
erosion modulus in the study area was 54.19 t/(hm-a) in 2005 and 37.96 t/(hm?.a) in 2015, which 
dropped by 29.95%. From the spatial distribution maps of soil erosion in 2005 and 2015, we 
found that the distribution patterns of soil erosion in both years were consistent. The areas of 
severe soil erosion were mainly concentrated in the south-central and northeast parts, in which 
soil erosion modulus was higher than 300 t/(hm?.a). In some of the areas with steep slopes, soil 
erosion modulus even exceeded 400 t/(hm?.a). Those areas are mostly deep canyon with broken 
land surface, high mountains, and deep valleys in the lower reaches of the Yangtze River and 
Yellow River. Soil erosion was relatively slight in the western and eastern parts with soil erosion 
modulus less than 100 t/(hm?-a), which was mainly due to the flat terrain and relatively high 
vegetation coverage. From the perspective of the severity of soil erosion, the areas with soil 
erosion modulus less than 100 t/(hm?.a) accounted for the most proportion (more than 80%) of the 
total area both in 2005 and 2015, followed by the areas of soil erosion modulus in the range of 
100-200 t/(hm?.a). Overall, the severity of soil erosion was moderate to high in the TRHR. 

From the aspect of dynamic changes in soil erosion, the total area of soil erosion increased 
from 12.01x10* km? in 2005 to 12.13x10*km? in 2015, but the average value of soil erosion 
modulus decreased. It can be seen from Figure 5d that, the reduction of soil erosion modulus was 
mainly concentrated in the range of 0-100 t/(hm?.a) in the study area, covering the largest 
proportion of the area with decreased soil erosion (308.16 km?). The area with the reduction of 
soil erosion modulus ranging from 100 to 200 t/(hm?.a) was the second largest (10.63 km?), which 
was mainly distributed in the middle south and northeast of the study area. The increased soil 
erosion modulus between 0 and 100 t/(hm?-a) covered the largest proportion of the area (8.48 
km?), followed by the range of over 200 t/(hm?.a). The areas with increased soil erosion modulus 
were mainly distributed in the southeast of the study area and had the characteristic of centralized 
distribution. The results illustrated that although the total area of soil erosion increased over the 
decade from 2005 to 2015, the overall severity of soil erosion was mitigated, presenting a trend of 
gradual improvement on soil erosion. 


3.3 Decomposition analysis of the C factor and the R factor influencing on the variation of 
soil erosion 


3.3.1 Impact of the C factor on soil erosion 


Vegetation is the most important environmental factor in controlling soil erosion, and it is also a 
positive factor in controlling soil erosion on slopes (Lieskovsky and Kenderessy, 2014). The 
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Fig. 5 Spatial distribution maps and statistical charts of soil erosion modulus in 2005 and 2015. (a), spatial 
distribution map of soil erosion modulus in 2005; (b), spatial distribution map of soil erosion modulus in 2015; 
(c), spatial distribution map of the variation of soil erosion modulus from 2005 to 2015; (d), area of different 
variation ranges of soil erosion modulus (the horizontal axis represents the different ranges of the variation in soil 
erosion modulus). 


influence of vegetation on soil erosion is mainly reflected in the prevention of surface runoff, 
which plays an important role in conserving soil. In this study, the C factor was determined by the 
NDVI, which is the most frequently used parameter in the study of the relationship between 
vegetation and soil erosion. We found that the average value of NDVI in the study area was 0.48 
and 0.44 in 2005 and 2015, respectively, indicating that vegetation in the study area had a slightly 
degrading trend. However, the vegetation coverage in some areas is increasing and the ecological 
conservation project is one of the main reasons for vegetation restoration (Liu et al., 2014). When 
the influence of the R factor is not considered, variation in soil erosion is only impacted by the C 
factor. Theoretically speaking, if vegetation coverage increases, soil erosion modulus will 
decrease correspondingly. 

The spatial distribution map of the variation of soil erosion caused by the C factor from 2005 to 
2015 was obtained using the additive decomposition model of the LMDI, and the results are 
shown in Figure 6. The positive and negative contribution values were categorized into three 
classes (0-100, 100-200 and >200 t/(hm?-a)) and the area of each class was calculated. As 
illustrated in Figure 6a and c, without considering other factors, the variation of soil erosion 
caused by the C factor was mainly distributed in the middle of the study area. The increased soil 
erosion modulus was mainly concentrated in the range of 0-100 t/(hm*-a) (covering an area of 
143.27x10° km?). The area where soil erosion modulus increased by more than 100 t/(hm?-a) due 
to the change in the C factor was relatively small. The area where soil erosion modulus decreased 
was mainly distributed in the west and east of the study area, and the decreased range was 
between 0 and 100 t/(hm?-a). We found that the decrease of soil erosion caused by the change in 
the C factor was greater than the increase of soil erosion during the period from 2005 to 2015, 
which suggested that vegetation can alleviate soil erosion in the TRHR. 

To understand the impact of the change in vegetation on soil erosion in different changing 
ranges of NDVI, we used the spatial analysis tool in the ArcGIS to overlay the NDVI changing 
map and the spatial variation map of soil erosion caused by the C factor. Then, we obtained the 
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cumulative contribution values of the C factor to the variation of soil erosion in different NDVI 
variation intervals (Fig. 7a). According to the results shown in Figure 7a, the cumulative 
contribution values of the C factor to the variation of soil erosion modulus had roughly the same 
variation characteristics in the increased and decreased ranges of NDVI. In the increased range of 
NDVI, the cumulative contribution values of the C factor to soil erosion was negative. With the 
increase of NDVI, the overall cumulative value showed an increasing trend first and then a 
decreasing trend. The cumulative contribution values of the C factor to the decreased soil erosion 
modulus was the highest when the increased NDVI was in the range of 0.05—-0.10, followed by the 
range of 0.10-0.15. Similarly, the cumulative contribution values of the C factor to the increased 
soil erosion modulus was the highest when the NDVI was in the range from —0.10 to —0.05. 
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Fig. 6 Spatial distribution maps (a and b) and statistics of the areas (c and d) of the contributions of the C factor 
and the R factor to the variation of soil erosion. (a), spatial distribution map of the contribution of the C factor to 
the variation of soil erosion; (b), spatial distribution map of the contribution of the R factor to the variation of soil 
erosion; (c), area chart of the contribution of the C factor to the variation of soil erosion; (d), area chart of the 
contribution of the R factor to the variation of soil erosion. 
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Fig. 7 Cumulative soil erosion contribution values of the C factor (a) and the R factor (b). (a), contribution 
values of the NDVI in different changing ranges to the variation of soil erosion from 2005 to 2015; (b), 
contribution values of the R factor in different changing ranges to the variation of soil erosion from 2005 to 2015. 
The vertical axes represent soil erosion contribution values of the NDVI change and the R factor change. 
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3.3.2 Impacts of the R factor on soil erosion 


The average annual precipitation in the TRHR was 476.70 mm in 2005 and 369.05 mm in 2015, 
and the average rainfall erosivity was 562.86 MJ-mm/(hm?-h-a) in 2005 and 381.86 
MJ-mm/(hm?-h-a) in 2015 (Figs. 2 and 4). The precipitation in 2015 was overall less than that in 
2005, and so was the rainfall erosivity. When the C factor was not considered, the variation of soil 
erosion was only impacted by the R factor. On the whole, the change in the R factor during the 
period from 2005 to 2015 mitigated soil erosion in the TRHR. According to the results of the 
LMDI model (Fig. 6b and d), soil erosion in some areas became more severe due to the change in 
the R factor. The contribution values of soil erosion caused by the R factor declined overall. The 
reduced soil erosion modulus caused by the R factor was mainly concentrated in the range of 0- 
100 t/(hm-a), covering an area of about 307.64 km?. The areas of which the reduction of soil 
erosion modulus between 100 and 200 t/(hm?-a) were mainly distributed in the central part. The 
area with increased soil erosion modulus caused by the R factor was mainly distributed in the 
southeast of the TRHR. Moreover, the area where soil erosion modulus increased was much 
smaller than the area where soil erosion modulus decreased (Fig. 6d). Precipitation in the study 
area is mainly concentrated in summer and autumn. The uneven spatial and temporal distribution 
of precipitation is one of the most important reasons for the misdistribution of the spatial 
variation of soil erosion. 

The spatial analysis tool in the ArcGIS was used to overlay the R factor change map and the 
spatial variation map of soil erosion caused by the R factor. Then, we obtained the cumulative 
contribution values of the R factor to the variation of soil erosion in different R variation intervals 
(Fig. 7b). The cumulative contribution values of the R factor to the variation of soil erosion in the 
decreased intervals of the R factor were much larger than those in the increased intervals. As the 
range increases, the cumulative contribution values showed an increasing trend first and then a 
decreasing trend. The comprehensive analysis reveals that the influence of the R factor on the 
variation of soil erosion in this area has a gradient-like property. 

3.3.3 Assessment of the dominant factors influencing soil erosion 


As is shown in Figure 8, soil erosion in the TRHR decreased greatly from 2005 to 2015. The 
reasons for the decrease in soil erosion modulus varied in different areas. Decreased soil erosion 
affected by the combination of the C factor and the R factor (DECR) occupied the largest 
proportion, accounting for 54.95% of the area where soil erosion decreased, and these areas were 
mainly distributed in the west and east of the study area. The impact of the R factor on the 
decrease of soil erosion was the second greatest (41.31% of the area), and the area was mainly 
distributed in the middle of the TRHR. The decreased soil erosion due to the C factor (DEC) 
covered a relatively small area (12.10x10? km), which was mainly distributed in the southeast of 
the study area. The rainfall erosivity was calculated from precipitation, reflecting the high 
correlation between rainfall erosivity and precipitation itself. To some degree, the relationship 
between soil erosion and precipitation can be understood as the correlation between soil erosion 
and the R factor. The area where soil erosion increased (12.95x10? km?) was smaller than the area 
where soil erosion decreased (325.01x10% km?) (Figs. 8 and 9). Among the influencing factors of 
increased soil erosion, the percentages of soil erosion caused by the C factor, the R factor, and the 
combination of the two accounted for 64.10%, 12.02%, and 23.88%, respectively. The C factor 
was the dominant factor causing the increase of soil erosion and the increased areas were mainly 
distributed in the middle and the east of the study area (Fig. 9). The areas where soil erosion 
increased due to the combined effects of the C factor and the R factor were mainly distributed in 
the south of the study area (Fig. 9). 


4 Discussion 


4.1 Advantages of the LMDI method in analyzing the influencing factors of soil erosion 


In this study, we carried out the research on soil erosion based on the LMDI method in the TRHR. 
There are several advantages using the LMDI method to examine the quantitative contribution of 
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Fig. 8 Spatial distribution map of dominant factors influencing the decrease of soil erosion. DEC, decreased soil 
erosion due to the C factor; DER, decreased soil erosion due to the R factor; DECR, decreased soil erosion due to 
the C factor and the R factor. 
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Fig. 9 Spatial distribution map of dominant factors influencing the increase of soil erosion. IEC, increased soil 
erosion due to the C factor; IER, increased soil erosion due to the R factor; IECR, increased soil erosion due to 
the C factor and the R factor. 


each influencing factor to soil erosion compared with the traditional methods, such as the 
statistical or regression approaches (Sun et al., 2013; Yao et al., 2016; Zhao et al., 2017; Zhang et 
al., 2018). First, the impact analysis of each influencing factor is independent of each other. As 
shown in Equations 8 and 9, there is no direct dependence between the C factor and the R factor. 
The analysis of each influencing factor was conducted by an independent formula, which does not 
affect the other factors. Second, there is no requirement for the order of factor decomposition in 
the calculation process. Third, the LMDI method can quantify the impacts of the influencing 
factors (Ang, 2015) and quantitatively calculate the contribution values of the influencing factors 
to the variation of soil erosion pixel by pixel, thus producing spatial distribution maps of the 
contribution values; in contrast, the regression relationships established by the traditional methods 
are mostly used to describe the influences of vegetation and precipitation on soil erosion in the 
whole region, which cannot quantitatively identify the specific contribution value from the pixel 
scale. Because the influences of vegetation and precipitation on soil erosion vary in different 
regions, the spatial heterogeneity would be neglected by directly establishing regression 
relationships in traditional methods. Besides, the traditional methods rely on experimental or 
statistical data that are difficult to obtain, especially in undeveloped areas. However, the LMDI 
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method can carry out the residual-free decomposition of the influencing factors of soil erosion 
based on the RUSLE model, using remote sensing data (freely available, such as MODIS and 
Landsat) and meteorological data instead of experimental or statistical data that are difficult to 
obtain. In this study, we have got the quantitative contribution values of the C factor and the R 
factor to the variation of soil erosion, and the results are shown in Figure 6. We can clearly 
understand the influences of vegetation and precipitation on soil erosion through the contribution 
maps of the C factor and the R factor. 


4.2 Factors influencing soil erosion in the TRHR 


Numerous studies have demonstrated that the variation of soil erosion is caused by the combined 
actions of vegetation and precipitation (Sun et al., 2013; Yin et al., 2014; Guo et al., 2016; He et 
al., 2018; Meng et al., 2018; Wang et al., 2018). The two factors have an interactive coupling 
relationship, whereas the interaction mechanism is still unclear. Grasping the influence 
mechanism of the influencing factors on the variation of soil erosion is the foundation of 
formulating a reasonable policy for protecting the ecological environment. The TRHR is located 
in the alpine region, where the ecosystem is extremely vulnerable (Han et al., 2018). In 2005, the 
government invested 7.5x10° CNY to launch the first stage of the ecological conservation and 
restoration project in the TRHR (Shao et al., 2017; Zhang et al., 2017), and the project had a 
positive impact on the vegetation recovery from 2005 to 2015 (Shen et al., 2018). As can be seen 
from Figure 5, soil erosion in the TRHR has been alleviated from 2005 to 2015 overall and the 
combination of the C factor and the R factor is the main reason for the decrease in soil erosion 
(Fig. 8), indicating the effectiveness of the project to some extent. The results are consistent with 
the previous studies (Shao et al., 2017; Shen et al., 2018). Ecological reconstruction and 
restoration efforts and climate change benefit ecosystem services. 

Decreased precipitation led to the decline in the R factor value from 2005 to 2015, which is a 
reason for the decrease in soil erosion (Li et al., 2014; Jiang et al., 2016). At the same time, the 
positive effects of precipitation on vegetation growth and restoration have been impaired (Zhang 
et al., 2018). Therefore, it is recommended to scientifically deploy artificial rain enhancement 
projects in the future (Jiang et al., 2016; Zhang et al., 2016; Shao et al., 2017). Moreover, 
precipitation is highly uncertain, which increases the difficulty of reducing soil erosion (Nearing 
et al., 2005). The IPCC have reported that the probability of precipitation may increase in the 
future (IPCC, 2018). Thus, severer events of soil erosion may occur in the TRHR due to the 
increased precipitation in the future, so it is necessary to pay more attention to soil erosion and 
take measures in advance. 

Vegetation has a significant effect on the variation of soil erosion, especially in the southeast 
corner of the study area (Figs. 8 and 9). Improving vegetation coverage is an effective way to 
reduce the degree of soil erosion (Zhang et al., 2016; Shen et al., 2018). However, due to the 
harsh ecological environment in the TRHR, the restoration of the root layer, which is of great 
significance for soil and water conservation, is extremely slow, and the recovery of soil physical 
and chemical properties is even slower (Shao et al., 2017; Xiong et al., 2019). This phenomenon 
also suggests that ecological restoration in the TRHR is a long-term and arduous task that requires 
continuous efforts (Huang et al., 2018). The main factors causing the change in vegetation are 
climate conditions and human activities (Zhang et al., 2016; Mishra and Mainali, 2017). Human 
activities such as overgrazing could cause vegetation to degrade, which would increase soil 
erosion accordingly. In general, ecological projects such as returning farmlands to forests and 
grasslands can reduce soil erosion by changing land surface coverage. Therefore, we suggest that 
the ecological conservation and restoration project should be continued to strengthen to reduce 
the damage caused by soil erosion. 

In the TRHR, the interaction of soil erosion and grassland degradation forms a vicious cycle of 
alpine grassland ecosystem degradation (Li et al., 2014). Moreover, rampant rodent infestation 
not only reduces the content of soil organic carbon, but also worsens the quality of soil (Shao et 
al., 2017). Grazing also has an important influence on soil properties (Fan et al., 2010; Lin et al., 
2017; Zhang et al., 2017). Soil-grass-animal interaction and mutual restriction make the 
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mechanism of soil erosion more complicated (Lin et al., 2017; Cao et al., 2018). In the future, 
measures need to be conducted to strengthen the scientific planning of grassland grazing and 
rodent pest management, prevent the aggravation of soil erosion caused by overgrazing and 
rodents, and maintain the balanced and stable development of the whole alpine ecosystem. 


4.3 Uncertainty analysis 


The mechanism of soil erosion in the TRHR is quite complicated, as there are various types of 
soil erosion including freeze-thaw erosion, wind erosion, and water erosion (Zhang et al., 2016; 
Shao et al., 2017). However, the RUSLE model mainly considers water erosion, which has 
become a technical defect of the popularization and application of the RUSLE model in a 
dynamic environment with various types of soil erosion (Lin et al., 2017). 

A series of methods are currently available for the spatial interpolation of precipitation data. In 
the process of calculation, the LMDI method relies on the spatial distribution map of 
precipitation, which is obtained from meteorological data through spatial interpolation. The 
interpolation accuracy of precipitation has a great influence on the results. In this study, spatial 
interpolation of precipitation was done by using ANUSPLIN, which has been extensively applied 
for spatially interpolating hydrometeorological surfaces all around the world (Hijmans et al., 
2005; McVicar et al., 2007; Zhang et al., 2013; Fick and Hijmans, 2017). As shown in Figure 2, 
we have got a relatively high interpolation accuracy of precipitation in both years (R?=0.90 in 
2005 and R?=0.95 in 2015). Besides, the number of meteorological stations in the Qinghai-Tibet 
Plateau is relatively small, exerting some impacts on the interpolation accuracy. Different 
interpolation methods might lead to different interpolation results, so future work would be 
focused on the issue. 

We have illustrated that the land use types remained virtually unchanged during the study 
period of 2005-2015, and the topography and soil properties may have some changes during this 
period. The variations in soil physical-chemical properties and soil texture and structure will lead 
to the change in the K factor (Sun et al., 2013; Jiang and Zhang, 2016), which can cause the 
variation in soil erosion. However, we did not have the DEM and soil data in different years due 
to the data limitation, so we could not analyze the influences of the LS factor and the K factor on 
the variation of soil erosion. 

More experimental or statistical data are required for a precise result. Due to the condition 
limitation, we did not conduct field experiments in the study area to verify the results of this 
study, e.g., the application of 1°’Cs and other nuclear tracer techniques in the ground verification. 
To validate the spatial distribution of the simulated soil erosion across the TRHR, we compared 
the estimated results by using the RUSLE model with the measured results or the results from 
previous researches (Huang et al., 2011; Shao et al., 2011; Teng et al., 2018), and the comparison 
showed that the results were consistent. Due to the lack of statistics of soil erosion in 2005, we 
compared soil erosion area (12.0110* km?) calculated by the RUSLE model with that (11.48x 104 
km?) from previous literature (Liu et al., 2014), and found that the overall results were also 
consistent. In spite of some uncertainties, the statistics from the government (Qinghai Provincial 
Water Resources Department; http://slt.ginghai.gov.cn/), together with some previous literature, 
can prove that the estimation of this study is satisfactory. 

The ecological conservation and restoration project in the TRHR was implemented in 2005 
(Zhang et al., 2017). By 2015, the first stage of the ecological project has been completed, and the 
second stage has already been started. The project aims to improve vegetation coverage and 
ecological quality by protecting grasslands, forests, and wetlands. A large number of sloping 
farmlands in the study area were artificially converted into grasslands and forests, and the 
vegetation coverage in local areas increased greatly (Shao et al., 2017), which significantly 
reduced the severity of soil erosion, reflecting the effects of the ecological conservation and 
restoration project in reducing soil erosion. However, there are still some areas showing the 
deterioration of vegetation. The ecological conservation and restoration project is a long process 
that needs to be continuously strengthened. 
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5 Conclusions 


In this study, the RUSLE model was used to evaluate soil erosion in the TRHR during the first 
stage of the ecological and restoration project between 2005 and 2015. To explore the influencing 
factors of soil erosion, we designed a method based on the LMDI model to quantitatively evaluate 
the contribution values of the C factor and the R factor to soil erosion. Soil erosion showed a 
downward trend overall, and the average soil erosion modulus decreased from 54.19 t/(hm?.a) in 
2005 to 37.96 t/(hm?.a) in 2015. The area where soil erosion amount increased was smaller than 
the area where soil erosion decreased. The cumulative contribution values of the C factor to the 
decrease of soil erosion were greater than those to the increase of soil erosion, which benefited 
from the improved vegetation coverage, indicating the effectiveness of the ecological 
conservation project. The dominant factor causing the decrease of soil erosion is the combination 
of the C factor and the R factor, followed by the R factor only; and the dominant factor 
contributing to the increase of soil erosion is the C factor. 

The method used in this study can provide a new idea for the further studies on exploring the 
influencing factors of soil erosion on a larger scale in the world, which has a broad prospect of 
popularization and application. A comprehensive assessment involving the LS factor and the K 
factor is needed to be carried out in the further studies. 
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Supplementary material 


1 Spatial distribution of total annual rainfall 
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Fig. S1 Spatial distribution of total annual rainfall in 2005 (a) and 2015 (b) 


2 Calculations of the RUSEL (Revised Universal Soil Equation) factors 


2.1 Vegetation cover factor (C factor) 


The C factor is closely related to vegetation coverage, which significantly influences soil erosion 
(Fu et al., 2011). We used the method developed by Cai et al. (2000) to calculate the C factor. 


1, f =0 
C =< 0.6508 — 0.34361gf, 0< f < 78.3%, (S1) 
0, f > 78.3% 
f - NDVI — NDVI y, (82) 


NDVI na, — NDVI pin 
where C is the vegetation cover factor (dimensionless); f is the vegetation coverage determined 


using Normalized Difference Vegetation Index (NDVI) derived from the SPOT-VGT NDVI; 
NDVImin is the NDVI of bare soil; and NDVImax is the regional maximum NDVI. 


2.2 Rainfall erosivity factor (R factor) 


The R factor reflects the detachment of soil at a location by raindrop splash erosion (Dissanayake 
et al., 2019). In this study, we used the following equation developed by Wischmeier and Smith 
(1978) to calculate the R factor. 


R= 31.735 = 1 Qt 5losio (P/P)-0.08188) (S3) 
i=1 
where R is the rainfall erosivity factor (MJ-mm/(hm?-h-a)); P; is the monthly rainfall (mm); and P 
is the annual rainfall (mm). 
2.3 Topographic factor (LS factor) 
The LS factor constitutes slope length and slope steepness (Ganasri and Ramesh, 2016). The 
Digital Elevation Model (DEM) was utilized to estimate the LS factor, while this factor is quite 
complicated to be calculated directly using GIS software such as ArcGIS. Considering the 
characteristics of soil erosion in China, we used the calculation tool of topographic factors 
developed by Fu et al. (2015), which is suitable for China and easy to use with a friendly 
interface. The main equations for the algorithm are shown as follows: 
m=0.2, 0<0.5° 
_( à \"Jm=0.3, 0.5°<@<1.5° 
m=0.4, 1.5°<8<3.0” 
m=0.5, 0 = 3.0° 


(S4) 
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10.8sin 0 + 0.03, 0<5° 
S= 16.8sin0—0.50, 5°<@<10°, (S5) 
21.91sin 6 — 0.96, 0 <10° 


where L is the slope length factor (dimensionless); S is the slope steepness factor (dimensionless); 
@ is the angle of the slope (°); m is a dimensionless constant; and à is the slope length (m), which 
depends on the steepness of the slope in percent. 

2.4 Soil erodibility factor (K factor) 


The K factor depends on the soil characteristics such as soil texture and organic matter content 
(Fayas et al., 2019). In this study, the K factor was determined using the EPIC (Erosion 
Productivity Impact Calculator) equation (Sharpley and Williams, 1990) as follows: 


0.3 
K =; 0.2 +0.3exp} —0.0256 x SAN x| 1 ate x SIL 
100 CLA +SIL 


s (S6) 


das 0.25 x OM das 0.7 xSN 
` OM + exp(3.72 — 2.950M) ` SN +exp(—5.51+ 22.9SN) 
where SAN, SIL, CLA, and OM are the percentage content of sand, silt, clay, and organic matter 
(%), respectively; and SN is equal to 1-SAN/100. 


2.5 Support practice factor (P factor) 


The P factor is the ratio of soil loss with a specific support practice to the corresponding loss with 
up slope and down slope cultivation (Dissanayake et al., 2019). The value of the P factor ranges 
from 0 to 1, in which O represents high-quality conservation practices and 1 indicates poor 
preservation conservation practices (Kang et al., 2018). 

Considering the land use types of the study area and the previous literature, the P factor value is 
0.15 for the agricultural land (Li and Zheng, 2012). There is no erosion in water (Fu and Zha, 
2008), wetland, bare rock, snow, and ice, thus the P factor value is 0 in these land use types. There 
are rarely water conservation measures in the rest land use types, so the P factor value is 1 (Zhao 
et al., 2007; Fu and Zha, 2008). 
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